use "$data/ctbdata_estimation", replace
* Estimate Individual-level parameters: 
* Regress the logged allocation ratio for the chosen option on a
* present period indicator t0
* the delay k in days
* the logged price ratio of the options in the budget
su id, det
forvalues i = 1(1)`r(max)' {
	di "Individual `i'"
	cap reg ln_consrat_star t0 k ln_pratio if id==`i' & choice_sd!=0 , noconstant 
	cap nlcom (a: 1/_b[ln_pratio] + 1) (b: exp(_b[t0]/_b[ln_pratio])) (d: exp(_b[k]/_b[ln_pratio])), post
	cap replace a_ind_ctb = _b[a] if id == `i' 
	cap replace b_ind_ctb = _b[b] if id == `i' 
	cap replace d_ind_ctb = _b[d] if id == `i' 
	cap replace a_ind_ctb_se = _se[a] if id == `i' 
	cap replace b_ind_ctb_se = _se[b] if id == `i'
	cap replace d_ind_ctb_se = _se[d] if id == `i' 
}

gen r_ind_ctb = (1/d_ind_ctb)^365 - 1
label var r_ind_ctb "Estimate: Annualized Discount Rate"

*** Aggregate-level parameter estimates using OLS
* Prepare variables: estimates and SEs to be populated

gen a_ctb=.
label var a_ctb "Estimate: Risk Aversion: Alpha"
gen b_ctb=.
label var b_ctb "Estimate: Present Bias: Beta"
gen d_ctb=.
label var d_ctb "Estimate: Long-run Discounting: Delta"
gen r_ctb=.
label var r_ctb "Estimate: Annualized Discount Rate"

gen a_ctb_se=.
label var a_ctb "Standar Error: Alpha"
gen b_ctb_se=.
label var a_ctb "Standar Error: Beta"
gen d_ctb_se=.
label var a_ctb "Standar Error: Delta"
gen r_ctb_se=.
label var a_ctb "Standar Error: Annualized Discount Rate"

* Loop over payoff treatment conditions
forvalues j=1/4 {
di "***** Payoff Treatment Condition=`j'****"
reg ln_consrat_star t0 k ln_pratio if ctb==`j', noconstant vce(cluster subject_id)

nlcom (d: exp(_b[k]/_b[ln_pratio])) (a: 1/_b[ln_pratio] + 1) (b: exp(_b[t0]/_b[ln_pratio])), post

replace a_ctb = _b[a] if ctb==`j'
replace b_ctb = _b[b] if ctb==`j' 
replace d_ctb = _b[d] if ctb==`j' 
replace a_ctb_se = _se[a] if ctb==`j' 
replace b_ctb_se = _se[b] if ctb==`j' 
replace d_ctb_se = _se[d] if ctb==`j' 

nlcom (r: 1/(_b[d]^365) - 1), post 
* Collect annualized discount rate
replace r_ctb =_b[r] if ctb==`j' 
replace r_ctb_se =_se[r] if ctb==`j'
}

* Results table aggregate estimates
*reg ln_consrat_star t0 k ln_pratio if ctb==1 , noconstant vce(cluster subject_id)
reg ln_consrat_star t0 k ln_pratio if ctb==1 & d_ind_ctb!=. , noconstant vce(cluster subject_id)
nlcom (d: exp(_b[k]/_b[ln_pratio])) (a: 1/_b[ln_pratio] + 1) (b: exp(_b[t0]/_b[ln_pratio])), post

* Compute implied discount factor based on annual rate reported in Andreoni et al 2015, Table 2 model 4
di (1/(1+.741))^(1/365)

* Annual rate
di (.998)^(-365)-1
* Discount rate
di (1/(1+1.076598))^(1/365)

* Keep one obs per respondent ans save to merge with other survey variables
bysort subject_id: gen aux=_n
keep if aux==1
save "$user/02 Data/ctb_ind.dta", replace

